Astronomy & Astrophysics manuscript no. high 


© ESO 2009 


March 26, 2009 





High-resolution Smoothed Particle Hydrodynamics simulations of 

the merger of binary white dwarfs 

p. Loren-Aguilar^'^, J. Isem^'^ and E. Garcfa-Berro^'^ 

' Departament de Ffsica Aplicada, Escola Politecnica Superior de Castelldefels, Universitat Politecnica de Catalunya, Avda. del 

Canal OKmpic 15, 08860 Castelldefels, Spain 
^ Institute for Space Studies of Catalonia, c/Gran Capita 2-4, Edif. Nexus 104, 08034 Barcelona, Spain 
' Institut de Ciencies de I'Espai, CSIC, Campus UAB, Facultat de Ciencies, Torre C-5, 08193 Bellaterra, Spain 

March 26, 2009 



ABSTRACT 



Context. The coalescence of two white dwarfs is the final outcome of a sizeable fraction of binary stellar systems. Moreover, this 
process has been proposed to explain several interesting astrophysical phenomena. 

Aims. We present the results of a set of high-resolution simulations of the merging process of two white dwarfs. 

Methods. We use an up-to-date Smoothed Particle Hydrodynamics code which incorporates very detailed input physics and an 

improved treatment of the artificial viscosity. Our simulations have been done using a large number of particles (~ 4 x 10') and 

cover the full range of masses and chemical compositions of the coalescing white dwarfs. We also compare the time evolution of 

the system during the first phases of the coalescence with that obtained using a simplified treatment of mass transfer, we discuss in 

detail the characteristics of the final configuration, we assess the possible observational signatures of the merger, like the associated 

gravitational waveforms and the fallback X-ray fiares, and we study the long-term evolution of the coalescence. 

Results. The mass transfer rates obtained during the first phases of the merger episode are in agreement the theoretical expectations. 

In all the cases studied the merged configuration is a central compact object surrounded by a self-gravitating keplerian disk, except in 

the case in which two equal-mass white dwarfs coalesce. 

Conclusions. We find that the overall evolution the system and the main characteristics of the of the final object are in agreement 
with other previous studies in which smaller resolutions were used. We also find the the fallback X-ray luminosities are of the order 
of 10*' erg/s. The gravitational waveforms are characterized by the sudden dissapearance of the signal in a few orbital periods. 

Key words, stars: white dwarfs — stars: interiors — stars: binaries: close — hydrodynamics — methods: numerical — accretion: 
accretion disks 



1. Introduction 

The coalescence of two close white dwarfs is thought to be 
one of the most common endpoints of the evolution of bi- 
nary systems. Consequently, the study of the coalescence pro- 
cess is an interesting issue, with many potential applications. 
Although the astrophysical scenarios in which the coalescence 
of two white dwarfs in a close binary system can occur and 
their relative frequencies have been well studied — see, for in- 
stance, Yungelson et al. (1994), Nelemans et al. (2001a, 2001b), 
and the recent review of Postnov & Yungelson (2006) — the 
merging process has received little attention until recently. The 
pioneering works of Mochkovitch & Livio (1989, 1990) who 
used an approximate method — the so-called Self-Consistent- 
Field method (Clement 1974) — and the full Smoothed Particle 
Hydrodynamic (SPH) simulations of Benz, Thielemann & Hills 
(1989), Benz, Cameron & Bowers (1989), Benz, Hills & 
Thielemann (1989), Benz et al. (1990), Rasio & Shapiro (1995) 
and Segretain, Chabrier & Mochkovitch (1997) were the only 
exceptions. Most of these early works had several drawbacks. 
For instance, some of them did not include a detailed nuclear 
network or the network was very simplistic, other used a very 
small number of SPH particles (~ 10^) and, finally, other did not 
discuss the properties of the merger configuration. Additionally, 
all these early works studied a reduced set of masses and chemi- 
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cal compositions and used the classical expression for the artifi- 
cial viscosity (Monaghan & Gingold 1983). This is an important 
issue since it is well known that SPH induces a large shear vis- 
cosity, which is more pronounced when the classical expression 
for the artificial viscosity is used. Additionally, energy dissipa- 
tion by artificial viscosity can lead to overheating and, hence, the 
peak temperatures achieved during the merger and the associated 
nucleosynthesis depend on the choice of the artificial viscos- 
ity. However, the situation has changed recently. Guerrero et al. 
(2004) opened the way to more reaUstic simulations. In particu- 
lar, in these calculations the standard prescription of Monaghan 
& Gingold (1983) for the artificial viscosity was used but the 
switch originally suggested by Balsara (1995) was employed to 
partially suppress the excess of dissipation. More recently, the 
simulations of Yoon et al. (2007) were carried out using a mod- 
em prescription for the artificial viscosity with time dependent 
parameters (Morris & Monaghan 1997) which guarantees that 
viscosity is essentially absent in those parts of the fluid in which 
it is not necessary, but there are other prescriptions which are 
also suitable. In particular the prescription of Monaghan (1997), 
which is based in Riemann- solvers also yields excellent results, 
and does not result in overheating. This is the prescription we 
use in the present work. 

It is also interesting to realize that the number of particles 
used in these kind of simulations has increased considerably in 
the last years, according to the available computing power. For 
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instance, Segretain et al. (1997) used ~ 6 x 10"^ SPH particles 
to simulate the coalescence of a 0.9 + 0.6 Mq system. Later, 
Guerrero et al. (2004) studied a considerable range of masses 
and chemical compositions of the merging white dwarfs em- 
ploying a sizeable number of particles (~ 4 x 10'*). This range of 
masses and chemical compositions included 6 runs in which sev- 
eral inital configurations were studied, involving helium, carbon- 
oxygen and oxygen-neon white dwarfs. Loren-Aguilar et al. 
(2005) simulated the coalescence of a system of two equal-mass 
carbon-oxygen white dwarfs of 0.6 Mq. However, in this work 
only one high-resolution simulation was done. More recently, 
Yoon et al. (2007) have studied in detail the coalescence of a bi- 
nary system composed of two white dwarfs of masses 0.6 and 
0.9 Mq using 2 x 10^ SPH particles. However, only one simula- 
tion was presented in this work. It is thus clear that a thorough 
parametric study in which several white dwarf masses and chem- 
ical compositions are explored using a large number of SPH par- 
ticles and a more elaborated treatment of the artificial viscosity 
remains to be done. 

Possible apphcations of these kind of simulations include 
the double degenerate scenario to account for Type la supernova 
outbursts (Webbink 1984; Iben & Tutukov 1984) and the forma- 
tion of magnetars (King, Pringle & Wickramasinghe 2001). Also 
three hot and massive white dwarfs members of the Galactic 
halo could be the result of the coalescence of a double white 
dwarf binary system (Schmidt et al. 1992; Segretain et al. 1997). 
Additionally, hydrogen-deficient carbon and R Corona Borealis 
stars (Izzard et al. 2007; Clayton et al. 2007) and extreme he- 
lium stars (Pandey et al. 2005) are thought to be the consequence 
of the merging of two white dwarfs. Finally, the large metal 
abundances found around some hydrogen-rich white dwarfs with 
dusty disks around them can be explained by the merger of a CO 
and a He white dwarf (Garcia-Berro et al. 2007). Last but not 
least, the phase previous to the coalescence of a double white 
dwarf close binary system has been shown to be a powerful 
source of gravitational waves which would be eventually de- 
tectable by LISA (Loren-Aguilar et al. 2005). 

Depending on the mass ratio of both stars and on the ini- 
tial conditions of the binary system the fate of double white 
dwarf binary systems is a merging process due to the loss of 
angular momentum through gravitational wave radiation. Stars 
orbit each other at decreasing orbital separations until the less 
massive one overfills its Roche lobe and mass transfer begins. 
According to the initial conditions mass transfer proceeds ei- 
ther in an stable or a dynamically unstable regime. The stability 
of mass transfer is an important issue. If the mass transfer pro- 
cess is stable, mass will flow at relatively low accretion rates 
and the whole merging process could last for several million 
years. On the contrary, if mass transfer proceds in an unstable 
way, the whole merging process finishes in a few minutes. The 
difi'erence between the two cases relies on the ability of the bi- 
nary system to return enough angular momentum back to the 
orbit. In fact, there are two competing processes. On the one 
hand, the donor star is supported by the pressure of degener- 
ate electrons and, hence, it will expand as it loses mass, thus 
enhancing the mass-transfer rate. On the other, if orbital angu- 
lar momentum is conserved the orbit will expand as the donor 
star loses mass thus reducing the mass-transfer rate. The pre- 
cise trade-off" between both physical processes determines the 
stability of mass transfer. Guerrero et al. (2004) found that all 
the systems merged in a few hundred seconds, corresponding to 
mass transfer rates of ~ 10"^ Mq s ^ Since the Eddington rate 
is of the order of 10"^ Mq yr"' the most massive white dwarf 
cannot incorporate the material of the disrupted secondary in 



such a short timescale and, thus, the secondary forms a hot at- 
mosphere and a heavy keplerian disk around the primary. This 
has been challenged by the simulations of Motl et al. (2002) and 
D'Souza et al. (2006). These authors used a grid-based three- 
dimensional finite-diff'erence Eulerian hydrodynamical code and 
found that when the stars are co-rotating mass transfer is stable. 
Nevertheless, it should be noted that these simulations were done 
using simplified physical inputs. For instance, they used a poly- 
tropic equation of state. More importantly, grid-based methods 
are known to poorly conserve angular momentum. In any case, 
it is clear that to assess the stability of mass transfer large spatial 
resolutions are required given the degenerate nature of the donor 
star, since once mass transfer begins the radius of the secondary 
increases very rapidly thereby increasing the mass-loss rate. 

In the present paper we study the coalescence of binary 
white dwarfs employing an enhanced spatial resolution (4 x 10^ 
SPH particles) and a formulation of the artificial viscosity which 
very much reduces the excess of shear. The number of particles 
used in our simulations is one order of magnitude larger than 
those used in our previous simulations (Guerrero et al. 2004) 
and within a factor of 2 of those used in modem simulations 
(Yoon et al. 2007; Loren-Aguilar et al. 2005). This allows us 
to resolve smaller scale lengths — by a factor (10^/10^)'^-' * 2 
— than in Guerrero et al. (2004). Moreover, this is done for a 
broad range of initial masses and chemical compositions of the 
coalescing white dwarfs, in contrast with most modem simula- 
tions in which only a single coalescence was studied in detail 
(Yoon et al. 2007). In particular we study the following cases: 
0.3 + 0.5 Mq, 0.4 + 0.8 M0, 0.6 -i- 0.6 M©, 0.6 + 0.8 and 
0.6 + 1.2Mo. Although we have computed a large number of 
mergers we will only discuss in detail the results of the merger 
of a 0.6 + 0.8 Mq binary system. The main results of the rest 
of the simulations are only given in tabular form, but can be 
provided upon request. Accordingly, we will devote most of the 
paper to compare the results of this simulation with those previ- 
ously available, and we postpone the discussion of the effects of 
the initial conditions to a forthcoming publication. The paper is 
organized as follows. In section 2 a brief description of our SPH 
code is given. In section 3 we describe our results and in section 
4 we compare them with those of other authors. In section 5 we 
discuss our simulations. Specifically, we pay special attention to 
discuss the stability of the mass transfer episode. This is done in 
section 5.1, whereas in sections 5.2 and 5.3 the possible observa- 
tional signatures arising from the merging process are studied. In 
particular, we consider the gravitational wave pattern of the sev- 
eral mergers studied here (section 5.2) and the X-ray emission 
that might be expected from the early phases of the disk evolu- 
tion (section 5.3), while in section 5.4 the long-term evolution of 
the merger is discussed. Finally, in section 6 we summarize our 
major findings, we elaborate on the possible implications of our 
work and we draw our conclusions. 

2. Input physics and method of calculation 

We follow the hydrodynamic evolution of the binary sys- 
tem using a Lagrangian particle numerical code, the so-caUed 
Smoothed Particle Hydrodynamics. This method was first pro- 
posed by Lucy (1977) and, independently, by Gingold & 
Monaghan (1977). The fact that the method is totally Lagrangian 
and does not require a grid makes it specially suitable for study- 
ing an intrinsicaUy three-dimensional problem Uke the coales- 
cence of two white dwarfs. We will not describe in detail the 
most basic equations of our numerical code, since this is a well- 
known technique. Instead, the reader is referred to Benz (1990) 
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Fig. 1. Temporal evolution of the density for the coalescence of the 0.6 + 0.8 Mq double white dwarf binary system. The positions 
of the particles have been projected onto the xy plane (left panels) and in the xz plane (right panels). The units of positions and 
densities are, respectively, 10^ cm and 10^ g/cm^^. Times are shown in the right upper corner of each panel. These figures have been 
done using the visualization tool SPLASH (Price 2007). [Color figure only available in the electronic version of the article]. 



where the basic numerical scheme for solving the hydrodynamic 
equations can be found, whereas a general introduction to the 
SPH method can be found in the excellent review of Monaghan 
(2005). However, and for the sake of completeness, we shortly 
describe the most relevant equations of our numerical code. 

We use the standard polynomic kernel of Monaghan & 
Lattanzio (1985). The gravitational forces are evaluated using an 
octree (Barnes & Hut 1986). Our SPH code uses a prescription 
for the artificial viscosity based in Riemann-solvers (Monaghan 
1997). Additionally, to suppress artificial viscosity forces in pure 
shear flows we also use the viscosity switch of Balsara (1995). 
In this way that the dissipative terms are essentially absent in 
most parts of the fluid and are only used where they are really 
necessary to resolve a shock, if present. Within this approach, 
the SPH equations for the momentum and energy conservation 
read respectively 



dt 
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where p,y - (pi + pj)/2 and Fij = (F/ + Fj)/2, and F is a 
positive definite function which depends only on |r| and and on 
the smoothing kernel h, used to express gradient of the kernel 
VWij - Fijrij. The signal velocity is taken as v^j^ = c,- + cj - 
Avij ■ Cij and the rest of the symbols have their usual meaning. 
We have found that a - 0.5 yields good results. 

We have found that it is sometimes advisable to use a 
different formulation of the equation of energy conservation. 
Accordingly, for each time step we compute the variation of the 



internal energy using Eq. (j2]) and simultaneously calculate the 
variation of temperature using: 



dTi 

dt j^iiC,)jPiPj 



Vij-ViW(\rij\,h)+q^ 



(3) 



where ^visc includes the contribution of viscous dissipation, 
which is computed in a way analogous to that of Eq. Q. For 
regions in which the temperatures are lower than 6 x lO^'^ K or 
the densities are smaller than 6 X 10^ g/cm^ Eq. ^ is adopted, 
whereas Eq. ([3]) is used in the rest of the fluid. Using this pre- 
scription we find that energy is best conserved. 

The equation of state adopted for the white dwarf is the sum 
of three components. The ions are treated as an ideal gas but tak- 
ing into account the Coulomb corrections (Segretain et al. 1994). 
We have also incorporated the pressure of photons, which turns 
out to be important when the temperature is high and the density 
is small, just when nuclear reactions become relevant. Finally the 
most important contribution is the pressure of degenerate elec- 
trons which is treated integrating the Fermi-Dirac integrals. The 
nuclear network adopted here (Benz, Thielemann & Hills 1989) 
incorporates 14 nuclei: He, C, O, Ne, Mg, Si, S, Ar, Ca, Ti, Cr, 
Fe, Ni and Zn. The reactions considered are captures of a parti- 
cles, and the associated back reactions, the fussion of two C nu- 
clei, and the reaction between C and O nuclei. All the rates are 
taken from Rauscher & Thielemann (2000). The screening fac- 
tors adopted in this work are those of Itoh et al. (1979). The nu- 
clear energy release is computed independently of the dynamical 
evolution with much smaller time-steps, assuming that the dy- 
namical variables do not change much during these time-steps. 
Finally, neutrino losses have been also included according to the 
formulation of Itoh et al. (1996) for the pair, photo, plasma, and 
bremsstrahlung neutrino processes. 

Regarding the integration method we use a predictor- 
corrector numerical scheme with variable time step (Serna, 
Alimi & Chieze 1996), which turns out to be quite accurate. 
Each particle is followed with individual time steps. With this 
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Fig. 2. Temporal evolution of the temperature (in units of 10^ K) for the coalescence of the same binary system shown in Fig.[l] The 
positions of the particles have been projected onto the xy plane (left panels) and in the xz plane (right panels). These figures have 
been done using the visualization tool SPLASH (Price 2007). [Color figure only available in the electronic version of the article]. 



Run 


MwD 






Mej 




T 

^ max 




H 




Af 


P 

^nuc 






a) 


0.3+0.5 


0.62 


0.18 


0.12 


10-3 


6.0 X 10** 


6.0 X 10** 


0.2 


6.1 X 10- 


-3 


300 


1 X lO''^ 


4x 10^' 


9 X 103* 


0.14 


0.4+0.8 


0.92 


0.28 


0.12 


10-3 


6.5 X 10'* 


6.0 X 10** 


0.2 


6.0 X 10 


-3 


166 


1 X 10"*^ 


5 X 10^^ 


4 X 103' 


0.23 


0.6+0.6 


1.10 


0.10 


0.50 


10-3 


6.3 X 10** 


6.2 X lO** 


0.07 


5.5 X 10- 


-3 


514 





3 X 10^^ 


1 X lO'" 


0.27 


0.6+0.8 


1.10 


0.30 


0.30 


10-3 


1.6 X 10' 


8.7 X 10** 


0.2 


5.0 X 10- 


-3 


164 


1 X lO-*' 


3 X 1028 


6x 10^" 


0.33 


0.6+1.2 


1.50 


0.30 


0.30 


10-3 


1.0 X 10'° 


1.0 X 10' 


0.2 


4.4 X 10- 


-3 


122 


2 X lO'*" 


8 X 1036 


5 X 10"° 


0.52 



Table 1. Summary of hydrodynamical results. Masses and radii are in solar units, times in seconds and energies in ergs. The 
maximum temperature achieved during each simulation, Tpeak, and the temperature of the hot corona formed around the primary at 
the end of the simulations, r^ax, are discussed in the text. The rotational velocity is expressed in s-'. 



procedure the energy and angular momentum of the system are 
conserved to a good accuracy. To avoid numerical artifacts, we 
only use equal-mass SPH particles, as it was done in Yoon et al. 
(2007). This was not the case of the simulations of Guerrero et al. 
(2004) in which the masses of the SPH particles were different 
for each one of the coalescing white dwarfs. In order to achieve 
an equilibrium initial configuration we relaxed each individual 
model star separately, so the two coalescing white dwarfs are 
spherically symmetric at the beginning of our simulations, as it 
was the case in in all previous simulations of this kind. In all 
the cases initially the two white dwarfs are in a circular orbit at 
a distance larger than the corresponding Roche lobe radius of 
the less massive component. The systems are not synchronized 
because, at least in the stage previous to the coalescence itself, 
the time scale for loss of angular momentum due to the emis- 
sion of gravitational radiation is so small that it remains quite 
unlikely that there exists any dissipation mechanism able to en- 
sure synchronization (Segretain et al. 1997). This is the same ini- 
tial configuration adopted by Yoon et al. (2007), Guerrero et al. 
(2004) and Segretain et al. (1997). However, as already stated, 
we will study synchronized systems in a forthcoming publica- 
tion. To this system we add a very small artificial acceleration 
term which decreases the separation of both components. Once 
the secondary fills its Roche lobe this acceleration term is su- 
pressed. This procedure is quite similar to that adopted in all 



previous works (Guerrero et al. 2004; Yoon et al. 2007; Segretain 
et al. (1997). We adopt this instant as our time origin. 

The chemical compositions of the coalescing white dwarfs 
depend on the mass of each star White dwarfs with masses 
smaller than 0.45 Mq have pure He cores. For white dwarfs with 
masses within this value and 1 . 1 Mq we adopt the corresponding 
chemical composition, namely, carbon and oxygen, with mass 
fractions Xq - 0.4 and Xq = 0.6 uniformingly distributed 
through out the core. Finally, white dwarfs more massive than 
1.1 Mq have ONe cores of the appropriate composition (Ritossa 
et al. 1996). This is essential in studying the resulting chemical 
composition of the merger, as it wiU be shown in Sect. 3, and 
it is a clear improvement with respect to recent high-resolution 
simulations (Yoon et al. 2007) in which only a single system of 
two carbon-oxygen white dwarfs was studied. 

3. Results 

Fig.[T]shows the temporal evolution of the logarithm of the den- 
sity for the coalescence of the 0.6 + 0.8 Mq white dwarf binary 
system. In the left panels the positions of the SPH particles have 
been projected onto the equatorial plane and in the right panels 
onto the polar plane. Time (in seconds) is shown on the right up- 
per corner of each panel. As can be seen in the uppermost left 
panels, the initial configurations of both white dwarfs are sym- 
metric. Soon after, the less massive white dwarf fills its Roche 
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Disk Corona 



Run 


0.3+0.5 


0.4+0.8 


0.6+0.6 


0.6+0.8 


0.6+1.2 


0.3+0.5 


0.4+0.8 


0.6+0.6 


0.6+0.8 


0.6+1.2 


He 


0.94 


0.98 








5 X 10-5 


0.72 


0.68 








5 X 10-5 


C 


2.4 X 10"^ 


7 X 10-3 


0.4 


0.4 


0.39 


0.11 


0.12 


0.4 


0.4 


0.26 





3.6 X 10"^ 


9 X 10-3 


0.6 


0.6 


0.60 


0.17 


0.18 


0.6 


0.6 


0.66 


Ne 


9 X 10"" 


6 X lO-'o 








3 X 10-3 


3 X 10-13 


7 X 10-' 








0.07 


Mg 


4x 10-'"* 


5 X 10-11 








7 X lO-'' 


1 X 10-12 


1 X 10-5 








8 X 10-5 


Si 


1 X 10"" 


3 X 10-1* 








1 X 10-5 


1 X 10-11 


1 X 10-* 








5 X 10-5 


S 


2 X 10-23 


1 X 10-1" 








3 X 10-' 


1 X 10-10 


2 X 10-* 








5 X 10-5 


Ar 


< 10-30 


7 X 10-2' 








1 X 10-' 


5 X 10-' 


9 X 10-* 








4 X 10-5 


Ca 


< 10-30 


< 10-30 








8 X 10-' 


1 X 10-** 


5 X 10-* 








1 X 10-* 


Ti 


< 10-30 


< 10-30 








7 X 10-' 


2 X 10-1 


1 X 10-2 








1 X 10-* 


Cr 


< 10-30 


< 10-30 








8 X 10-' 


4x 10-* 


2 X 10-3 








2 X 10-* 


Fe 


< 10-30 


< 10-30 








5 X 10-* 


2 x 10-5 


1 X 10-5 








6 X 10-* 


Ni 


< 10-30 


< 10-30 








6 X 10-1 


2 X 10-' 


4 X 10-** 








1 X 10-2 


Zn 


< 10-30 


< 10-30 








6 X 10-* 


2x 10-' 


6x 10-1" 








2 X 10-5 



Table 2. Averaged chemical composition (mass fractions) of the heavy rotationally-supported disk and the hot corona obtained by 
the end of the coalescing process. 



lobe and mass tranfer begins, as can be seen in top central panel 
of this figure. The top right panel of Fig.[T]shows that, after some 
time, the matter outflowing the secondary hits the surface of the 
primary white dwarf and spreads on top of it. Note as well that 
since the radius of white dwarfs scales as ~ M-''-', as the sec- 
ondary loses mass its radius increases and, hence, the mass-loss 
rate of the secondary increases, thus leading to a positive feed- 
back of the process. As a consequence of this positive feedback 
an accretion arm is formed which extends from the remnant of 
the secondary white dwarf (central panels in Fig. [TJ to the sur- 
face of the primary white dwarf. This accretion arm becomes 
entangled as a consequence of the orbital motion of the coalesc- 
ing white dwarfs and adopts a spiral shape (bottom left panel). 
Ultimately, the secondary is totally disrupted and a heavy disk is 
formed around the primary (bottom central panel of Fig.[T]i. The 
bottom right panel of Fig.[T]shows that at time t - 152 s the disk 
is still not well formed and the remnant of a spiral arm still per- 
sists. We followed the evolution of this merger for some more 
time and we found that the final configuration has cylindrical 
symmetry, that most of the orbits of the SPH particles belong- 
ing to the secondary have been circularized and that the spiral 
pattern has totally disappeared. At the end of the simulations the 
radial extension of the disk is ~ 0.2 Rq, whereas its height is 
~ 5.0 X IQ-^Rq. 

The temporal evolution of the temperature for the merger of 
a 0.6 + 0.8 Mq binary system is shown in Fig. ^ As can be 
seen in this figure the material of the secondary is first heated 
by tidal torques. As the secondary begins the disruption process 
this material is transferred to the surface of the primary and, con- 
sequently, it is compressed and its temperature increases. The 
peak temperatures (Tpeak) achieved during the coalescence are 
displayed in Table [T] for each one of the runs presented in this 
paper. For the 0.6 + O.^Mq simulation the peak temperature 
is Tpeak ~ 1.6 X 10^ K, clearly larger than the carbon igni- 
tion temperature Tign ~ 10^ K, and occurs during the first and 
most violent part of the merger. However, a strong thermonu- 
clear flash does not develop because although the temperature in 
the region where the material of the secondary first hits the pri- 
mary increases very rapidly, degeneracy is rapidly lifted, lead- 
ing to an expansion of the material, which, in turn, quenches 
the thermonuclear flash. This is in agreement with the results of 
Guerrero et al. (2004) and Yoon et al. (2007). Thus, since these 



high temperatures are attained only during a very short time in- 
terval thermonuclear processing is very mild for this simulation. 
It is also interesting to compare the equatorial and polar distri- 
bution of temperatures shown in the central panels of Fig. |2] 
This comparison reveals that the heated material is rapidly re- 
distributed on the surface of the primary and, as a consequence, 
a hot corona forms around the primary. The spiral structure pre- 
viously described can be more easily appreciated in the bottom 
right panels of figure |2] In fact, this spiral structure persists for 
some more time. 

In all the cases studied here a self-gravitating structure forms 
after a few orbital periods, in agreement with our previous find- 
ings (Guerrero et al. 2004) and with those of Yoon et al. (2007). 
The time necessary for its formation depends on the system be- 
ing studied, and ranges from ~ 120 seconds to ~ 520 seconds. 
This self-gravitating structure consists in all the cases but that in 
which two equal-mass white dwarfs are involved in a compact 
central object, surrounded by a heavy keplerian disk of variable 
extension. In the case in which two 0.6 Mq white dwarfs are in- 
volved the configuration is rather different. In this last case the 
simmetry of the systems avoids the formation of a clear disk 
structure, giving rise instead to a rotating elipsoid around the 
central compact object, surrounded by a considerably smaller 
disk. In table [T] we summarize the most relevant parameters of 
all the mergers studied here. Columns two, three, four and five 
list, respectively, the mass of the central white dwarf obtained 
after the disruption of the secondary, the mass of the keplerian 
disk, the accreted and the ejected mass. All the masses are ex- 
pressed in solar units. In column six we show the peak tempera- 
tures achieved during the coalescence. In column seven we dis- 
play the temperature of the hot corona around the central object 
by the end of our simulations, whereas in column eight the ra- 
dius of the disk is shown. In column nine we list the disk half- 
thickness. Column ten displays the duration of the coalescence 
process and columns eleven, twelve and thirteen display the en- 
ergetics of the process. In particular we show the thermonuclear 
energy released during the coalescence process {E„^c), the neu- 
trino energy {Ey) and the energy radiated in the form of gravita- 
tional waves (Egw)- Finally, in column fourteen we list the an- 
gular velocities of the central compact remnants. As can be seen, 
for the first two simulations the accreted mass is approximately 
the same and the same occurs for the last two simulations. As al- 
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Fig. 3. Left panels: rotational velocity of the merger products as a function of the radius. For the sake of comparison the keplerian 
velocity is also shown as a dashed line. Right panels: surface density profiles compared with the theoretical thin disk model profiles 
(dashed lines). 



ready commented, the simulation in which two equal-mass white 
dwarfs are involved is rather special and in this case we do not 
have properly speaking a disk, although a very flattened region 
with cylindrical symmetry forms around a central object of el- 
lipsoidal shape. The mass of this region is ~ 0.1 Mq. In all five 
cases the mass ejected from the system (those particles which 
acquire velocities larger than the escape velocity) is very small 
(~ 10"^ Mq), and thus the merging process can be considered 
as conservative. The maximum temperatures of the coronae in- 
crease as the total mass of the binary system increases. It should 
be noted that for the case of the 0.6 -i- 0.6 Mq binary system 
the maximum temperature occurs at the center of the merged 
configuration. We have found that these temperatures are some- 
what smaller than those obtained in our previous simulations 
(Guerrero et al. 2004). This is a direct consequence of the im- 
proved treatment of the artificial viscosity and of the enhanced 
spatial resolution. For instance, we have found that using an en- 
hanced resolution and an improved prescription for the viscosity 
the peak temperature obtained in the simulation in which a 0.6 
and 0.8 Mq are involved is Tpeak - 1-64 x 10^ K — see Table 
1. When a reduced number of SPH particles and the classical 
expression for the artificial viscosity are used this temperature 
is T'peak - 1.72 X 10^ K, whereas when a reduced number of 
particles and an improved artificial viscosity are used the peak 
temperature turns out to be Tpeak - 1-68 x 10^ K. It is worth 
noting that the radial extension of the disks is roughly the same 
for all but one the simulations presented here and it is consider- 
ably smaller for the case in which two equal-mass white dwarfs 
are involved. This is a natural behavior since in this last case 
the central object is rather massive. Finally, it is as well interest- 
ing to realize that all the disks are rather thin, being the typical 
half-thickness of the order of ~ 10"^ Rq, much smaller than the 
typical disk radial extension, ~ 0.2 Rq. 

The chemical composition of the disk formed by the dis- 
rupted secondary can be found for all the simulations presented 
in this paper in Table |2] In this table we show, for each of 



the mergers computed here, the averaged chemical composition 
(mass fractions) of the heavily rotationally-supported disk — left 
section of table [2] — and the hot corona — right section — de- 
scribed previously. For the mergers in which two carbon-oxygen 
white dwarfs are involved the disk is mainly formed by carbon 
and oxygen and the nuclear processing is very small (see the 
peak temperatures shown in column ten of Table [T}. This is not 
the case for the simulations in which a lighter He white dwarf 
is involved. Since in these cases the Coulomb barrier is consid- 
erably smaller, the shocked material is nuclearly processed and 
heavy isotopes form. This is more evident for the case in which a 
massive He white dwarf of 0.4 Mq is disrupted by a massive CO 
white dwarf of 0.8 Mq — third and eight columns in Table |2] In 
this case the abundances in the disk and the hot corona are rather 
large. Note as well that the abundances of heavy nuclei in the 
hot corona are much larger than those of the disk, indicating that 
most of the nuclear reactions occur when the accretion stream 
hits the surface of the primary. Nevertheless, although the disk is 
primarily made of the He coming from the disrupted secondary 
the abundances of C and O are sizeable and, moreover, the disk is 
contaminated by heavy metals. This has important consequences 
because it is thought that some of the recently discovered metal- 
rich DA white dwarfs with dusty disks around them — also 
known as DAZd white dwarfs — could be formed by accretion 
of a minor planet. The origin of such minor planets still remains 
a mistery, since asteroids sufficiently close to the white dwarf 
would have not survived the AGB phase (Villaver & Livio 2007). 
However, planet formation in these metal-rich disks is expected 
to be rather efficient, thus providing a natural environment where 
minor planetary bodies could be formed and, ultimately, tidally 
disrupted to produce the observed abundance pattern in these 
white dwarfs (Garcfa-Berro et al. 2007). Nuclear reactions are 
also important in the case in which a regular 0.6 Mq carbon- 
oxygen white dwarf and a massive oxygen-neon white dwarf of 
1.2 M0 are involved. In this case the peak temperature achieved 
during the coalescence is rather high Tpeak - 1.0 x lO"' — see 
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Table [T] — enough to power carbon burning. Consequently, the 
chemical abundances of the keplerian disk and of the hot corona 
are largely enhanced in oxygen and neon, which are the main 
products of carbon burning. We must add, however, a caution- 
ary remark regarding the chemical compositions of the merg- 
ers studied here. White dwarfs are characterized in ~ 80% of 
the cases by a thin hydrogen atmosphere of ~ 10^* Mq on top 
of a helium buffer of ~ 10"^ Mq. In the remaining ~ 20% of 
the cases the hydrogen atmosphere is absent. Small amounts of 
helium could indeed change the nucleosynthetic pattern of the 
hot corona in all the cases studied here. Studying this possibility 
is beyond the scope of this paper and, thus, the changes in the 
abundances associated to burning of the helium buffer and of the 
atmospheric hydrogen remain to be explored. 

In figure|3]we explore the final characteristics of the merged 
configuration. We start discussing the left panels of figure [3] 
which show the rotational velocity of the merger as a function 
of the distance to the center of the merged object. Clearly in all 
the cases there is a central region which rotates as rigid solid — 
see last column of table [T] This behavior was already found in 
Guerrero et al. (2004) and Yoon et al. (2007), and it is a con- 
sequence of the conservation of angular momentum. On top of 
this region a differentially rotating layer is present. This rapidly 
rotating region is formed by material coming from the disrupted 
secondary, which has been accumulated on top of the primary 
and thus carries the original angular moment of the secondary. 
Finally, for sufficiently large radius a rotationally-supported disk 
is found. The exact location where the disk begins can be easily 
found by looking at the left panels of Fig.|3] where the keplerian 
velocity is also shown as a dashed line. The change in the slope 
of the profile of the rotational velocity clearly marks the outer 
edge of the compact inner object and the beginnig of the disk. 
All the disks extend up to some solar radii — see column eight 
in Table [T] The stratification of surface densities of these disks 
can be seen in the left panels of Fig. |3] where we have plotted 
the surface density as a function of the distance. For the sake 
of comparison the theoretical surface density of a thin disk an- 
alytical model (Livio, Pringle & Wood 2005) is also shown as 
a dashed line. Within this model the surface density of the disk 
should be of the form Z oc R^^. We have used /3 - 7/4 to pro- 
duce the dashed lines in the right panels of Fig. |3] very close to 
the value adopted by Livio et al. (2005), /3 - 3/2. As can be seen 
in this figure for the first two fiducial mergers studied here there 
is a region where the analytical model and the numerical results 
are in good agreement. However, at large enough distances the 
SPH density profile falls off more rapidly than that of the theoret- 
ical model. In the case of the merger of two equal-mass 0.6 Mq 
white dwarfs the agreement is poor In this case, the symmetry of 
the system avoids the formation of a clear disk structure, giving 
rise instead to a rotating ellipsoid around the central compact 
object. Moreover, it can be shown that the angular momentum 
of the disk can be expressed in terms of the disk radius R^isk 
and the disk mass Mjisk as 7. = ^Mdisk(GMwD^disk)'^^, where 
^ - (2 - /3)/(5/2 - P) - 1/3. The theoretical angular moments 
obtained using this equation agree very well with the results of 
our SPH simulations. 

In Fig. |4] we show the temperature profiles at the end of the 
simulations for some of the mergers studied here. We have av- 
eraged the temperatures of those particles close to the orbital 
plane. The average was done using cylindrical shells and the size 
of these shells was chosen in such a way that each of them con- 
tain a significative number of particles. As can be seen, for the 
0.3 -I- 0.5 Mq and the 0.4 + 0.8 Mq systems, the region of max- 
imum temperatures occurs off'-center, at the edge of the original 
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Fig. 4. Radially averaged temperature profiles as a function of 
radius. 



primary, in the region of accreted and shocked material, whereas 
for the merger in which two equal-mass 0.6 Mq white dwarfs 
coalesce the maximum temperature occurs at the center of the 
merged object, as it should be expected. These maximum tem- 
peratures are listed in the seventh column of Table [T] In fact, 
the temperature profiles shown in this figure clearly show that 
the cores of the primaries in the first two simulations remain al- 
most intact and, hence, are rather cold. These cores, in turn, are 
surrounded by a hot envelope wich corresponds to the shocked 
material coming from the disrupted secondary. Nuclear reactions 
are responsible for the observed heating of the accreted matter, 
initially triggered in the shocked regions. The case in which two 
0.6 Mq white dwarfs coalesce is somewhat different. In this case 
there is not a hot envelope around a central object — although a 
local maximum of temperatures can indeed be appreciated at the 
edge of the rapidly spinning central object, as shown in the bot- 
tom panel of Fig.|4] — and, instead, the central region of the com- 
pact object is formed by the cores of the merging white dwarfs. 
Most of the temperature increase in this case is due to viscous 
heating since nuclear reactions are negligible because the in- 
crease in temperature of the shocked material is not enough to 
ignite carbon. In all the cases it is rather apparent a sizeable dis- 
persion of temperatures in the outermost regions. This dispersion 
is due in part to the fact that this region contains some of the par- 
ticles that were ejected during the first and most violent phases 
of the merging process. 

In summary, we have found that the qualitative behavior of 
all the mergers studied here is similar. In particular the less mas- 
sive component of the binary system is disrupted in a few orbital 
periods. Additionally, we find that nuclear reactions are only im- 
portant in the cases in which the secondary is a He white dwarf 
or in the case in which the primary (the accretor) is an ONe 
white dwarf. The total ejected mass is very small in all the cases. 
Finally, the overall final configuration is very similar in all the 
cases but in the case in which two equal-mass white dwarfs coa- 
lesce. 
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4. Comparison with previous worlds 



As previously discussed, all the mergers studied here coalesce in 
a dynamical time scale, regardless of its mass ratio. Additionally, 
in agreement with previous calculations, the central regions of 
the remnant rotate as a rigid body. On top of this rapidly spin- 
ning core a hot corona forms and, on top of it, a heavy keple- 
rian disk. The question is now, how the main characteristics of 
the merged configuration compare with those obtained in pre- 
vious works? Yoon et al. computed the coalescence of a single 
0.6 -I- 0.9 Mq binary system and the comparison is not straight- 
forward, but our 0.6-1-0.8 Mq run is rather similar. The first thing 
to be noted is that the duration of the merger is very similar in 
both cases. We obtain a duration of 164 s and Yoon et al. (2007) 
obtain ~ 150 s. The central density of the rapidly spinning core 
is in both cases p ~ 10^ g/cm-*. The temperature of the core is 
log T = 7.9 whereas Yoon et al. (2007) obtain log T ~ 4.5, but 
this is due to our choice of the initial temperature of the coalesc- 
ing white dwarfs, for which we adopted T - 10^ K. However, the 
temperatures of the hot coronae are remarkably similar in both 
cases, log T ~ 8.5 and 8.6 respectively. The peak temperatures 
attained during the merger are also very similar — 1.7 x 10^ K 
and 1 .6x10'^ K, respectively. However, the temperature of the hot 
coronae is considerably smaller in their case T^ax - 5.8 x 10^ 
K. This value has to be compared with that shown in table [T] 
TmsLx - 8.7 X 10'^ K. However, it should be taken into account 
that Yoon et al. (2007) followed the evolution of the merger for 
much longer times. The sizes of the resulting disk are also very 
similar. Yoon et al. (2007) obtained 8 x 10"^ cm, whereas we 
obtain 1 x 10^° cm. Moreover, despite of the very different ap- 
proaches for the artificial viscosities adopted in the work of Yoon 
et al. (2007) and in the present work, the rotational velocities of 
the central spinning object are very close, 0.21 s ' in the case of 
Yoon et al. (2007), whereas we obtain 0.26 s"^ However, we find 
that the rotational velocities of the hot coronae are somewhat dif- 
ferent. In particular, Yoon et al. (2007) obtained acj ~ 0.54 s"' 
while we obtain 0.33 s ' . This could be due to the different treat- 
ment of the artificial viscosity and to the different masses of the 
coalescing white dwarfs. 



The comparison with the results of Guerrero et al. (2004) 
also yields interesting results. For instance, the angular veloc- 
ity of the central compact object of the 0.6 + 0.8 Mq merger in 
the calculations of Guerrero et al. (2004) is 0.33 s"', somewhat 
larger than that obtained here. Consequently, the central density 
obtained in the simulations of Guerrero et al. (2004) is smaller 
(p ~ 6.3 X 10^ g/cm-'). This is due to our improved treatment of 
the artificial viscosity which considerably reduces the excess of 
shear and, thus, translates into smaller centrifugal forces. Also, 
the position of the hot corona is smaller in the case of the simu- 
lations presented here. Specifically, the hot corona is located at 
0.008 /?0 in the simulations discussed here, whereas in the case 
of Guerrero et al. (2004) was located at 0.01 Rq. Also, the tem- 
perature of the central object is smaller in our case (T ~ 7.6x 10^ 
K) than in the case of Guerrero et al. (2004) for which a cen- 
tral temperature of T ~ 3.8 x 10^ K was obtained, even though 
the initial temperatures of the coalescing white dwarfs were the 
same in both cases (10^ K). Hence, our improved treatment of 
the artificial viscosity results in a smaller overheating and in 
smaller shear. 



5. Discussion 

5.1. Comparison with theory 

To obtain a better understanding of the coalescence process and 
to compare our results with those theoretically expected, we have 
numerically solved the equations of the evolution of the binary 
system during the mass transfer phase. The evolution of a bi- 
nary system during this phase is determined by three basic physi- 
cal processes, namely, gravitational wave emission, tidal torques 
and mass transfer There is a wealth of literature dealing with this 
problem. We adopt as our starting point the analysis of Marsh 
et al. (2004) and the more recent formulation of Gokhale et al. 
(2007). Within this approach the evolution of the orbital separa- 
tion a is given by 



a 7gw /i(Q-wi) /2(Q-W2) , M2 

TT- = -J -J -J ^^'^'I^JT 



(4) 



where Mi and M2 are, respectively, the masses of the accretor 
and of the donor, q = M2/M1 is the mass ratio, Q is the orbital 
velocity, iw, are the spin velocities, /, stands for the moments of 
inertia, r, are the synchronization timescales and q-^ is defined as 



q^=l+M2 



]2- ]\ 



(5) 



orb 



being j\ the specific angular momentum of the matter arriving to 
the accretor and j2 the angular momentum of the matter leaving 
the donor In the calculations presented here we have adopted for 
ji the expression for disk fed accretion: 



./I 



^JGMlRl 



whereas for 72 we have 



J2 



(6) 



(7) 



The first term in Eq. Q corresponds to the change in the orbital 
separation due to gravitational losses. Due to the short duration 
of the coalescing process, its contribution can be neglected. The 
second and third term describe the tidal couplings. Finally, the 
last term in Eq. (j4]) corresponds to the advected angular momen- 
tum. 

The evolution of the Roche lobe radius R]^ is given by 



2Rl 



Jew 

J orb 



/l(Q-Wi) /2(Q-W2) 



JovbTl 

2 



JmbT2 



M2 
M2 



(8) 



where ^r^^ ^ 0.30 + 0.16^, which is valid for 0.01 < ^ < 1 
(Gokhale et al. 2007). As can be seen from this set of equations, 
it s possible for a system to change the type of mass tranfer from 
stable to unstable or vice versa because q has a dynamical value. 
Thus, to determine the exact type of mass transfer of the system 
it is necessary to study how q evolves with time. 

To solve the previous equations, an expression for the mass- 
transfer rate is needed. The mass-transfer rate is essentially de- 
termined by the Roche lobe overfill factor, wich is defined as 
A = R2- Rl- For the Roche lobe radius we adopt the expression 
of Eggleton (1983). In the case of a 0.6 M0 donor white dwarf 
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a polytropic equation of state with n - 312 can be adopted and, 
thus, we have (Paczyriski & Sienkiewicz 1972): 



M2 = -aW(;u)^Aj 



being 



(9) 



(10) 



where ^ — M2l{M\ + M2) and a is a smoothly varying func- 
tion of the stellar parameters which, following the approach of 
Gokhale et al. (2007), we have taken as a constant freely ad- 
justable parameter. The synchronization timescales have been 
computed using the expressions of Campbell (1984) 



T2 



<M2 



(11) 



where, again, the normalization factors Tj and are freely ad- 
justable parameter (Marsh et al. 2004; Gokhale et al. 2007). With 
all these inputs the equations for the evolution of the binary sys- 
tem — Eqs. Q and ([8]l — together with the equations for the 
evolution of the spin angular velocities of each of the compo- 
nents 



(12) 



can be integrated. In doing so it has to be taken into account 
the logical limitations of the theoretical approach. In particular, 
the SPH results show that mass transfer is not perfectly con- 
servative, although this assumption is fairly good — see Table 
[T] Moreover, stars are not point-like masses and, most impor- 
tantly, we have not adopted in our analysis an equilibrium mass- 
radius relationship. These assumptions may produce marked dif- 
ferences between the SPH and the theoretical results. Perhaps the 
most critical assumption in determining the evolution of the sys- 
tem is the assumption that an equilibrium mass-radius relation 
holds for both members of the binary system. In fact, at the be- 
ginning of the mass-transfer episode stars are not in equilibrium. 
Consequently, we have adopted a different approach. In particu- 
lar when integrating Eqs. ( 12 1 we have computed for each time 
step the actual moment of inertia of each star In particular, in 
the case of the primary white dwarf for each computed model 
we look for location of the region with maximum temperature 
(see Fig.|2|. We then compute the mass interior to this shell and 
the corresponding moment of inertia. For the case of the donor 
white dwarf we look for the region which still has an approx- 
imate spherical symmetry (see Fig. [T]i and we follow the same 
procedure adopted for the accretor 

In figure 5]we show a comparison of the theoretical results — 
shown as a dashed line — and the SPH results — shown as dots 
— for the time evolution of the orbital separation and of the spin 
angular moment of the accretor (Ji) and donor stars (J2)- The 
three adjustable parameters adopted in the theoretical calcula- 
tions are, respectively, a = 14-5 Mq Rq /yr, = 3.50 x 10"* yr 

and = 4.7510"*^ yr As can be seen, during the first phases of 
the merger the agreement is excellent. Note however that we can 




Fig. 5. A comparison of the SPH values for the orbital distance a 
and for the spin angular moments of the donor and the accretor 
stars for the 0.6 + 0.8 Mq case and those obtained using a sim- 
plified model. The sping angular moments are expressed in units 
of 10^' g cm-/s. 



only compare the SPH results with the theoretical expectations 
while the secondary still preserves partially its initial shape. This 
is why we only show in Fig. |5] a reduced time interval, corre- 
sponding to the first five panels in Figs.[T]and|2] For times longer 
than ~ 70 s, the secondary rapidly dissolves and, hence, the ap- 
proach followed here is no longer valid. It is worth realizing that 
T*!* » t". This means that the synchronization timescale of the 
primary is much larger than that of the secondary. Accordingly, 
during this phase of the mass-transfer episode the donor rapidly 
synchronizes whereas the primary does not. Consequently, or- 
bital angular moment is transferred from the orbit to the donor 
in a short timescale, thus reducing the orbital separation. This, 
in turn, increases the mass-transfer rate and the final result is 
that the secondary is rapidly disrupted. Since the total angular 
momentum is conserved the material transferred to the primary 
must rotate rapidly, thus producing the characteristic rotational 
profiles shown in the left panels of Fig.|3] In summary, the results 
of the hydrodynamic calculations can be accurately reproduced 
by a simple model once all the weaknesses of the theoretical ap- 
proach are correctly taken into account. 

5.2. Gravitational wave radiation 

Gravitational wave radiation from Galactic close white dwarf 
binary systems is expected to be the dominant contribution to 
the background noise in the low frequency region, which ranges 
from ~ 10"^ up to ~ 10"^ Hz (Bender et al. 1998). Moreover, 
since during the merging process a sizeable amount of mass is 
transferred from the donor star to the primary at considerable 
speeds, the gravitational wave signal is expected to be detectable 
by LISA (Guerrero et al. 2004; Loren-Aguilar et al. 2005). It is 
thus important to characterize which would be the gravitational 
wave emission of the white dwarf mergers studied here and to 
assess the feasibility of dectecting them. 
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Fig. 6. Gravitational wave emission from the merger of a 0.6 + 
0.8 Mq close white dwarf binary system. The dimensionless 
strains h+ and are measured in units of 10"^'. The source 
has been assumed to be at a distance of 10 kpc. 



To compute the gravitational wave pattern we proceed as in 
Loren-Aguilar et al. (2005). In particular, we use the weak-field 
quadrupole approximation (Misner et al. 1973): 



2G d^QY(t-R) 



'jk_ 



(13) 



where t-R - t-d/c is the retarded time, d is the distance to the 
observer, and Qj^(t - R) is the quadrupole moment of the mass 
distribution, wich is given by 



pix, t - R)(x-'x'' - ^x-'Sjk)d''x 



(14) 



To calculate the quadrupole moment of the mass distribution us- 
ing SPH particles, Eq. ( 14 1 must be discretized according to the 
following expression 



Using this prescription, the corresponding strains for the 
0.6 + 0.8 Mq, which is a representative case, are shown in Fig.|6] 
As shown in this figure the gravitational waveforms rapidly van- 
ish in a couple of orbital periods and the gravitational wave emis- 
sion during the coalescence phase does not have a noticeable 
large peak. Hence, the gravitational wave emission is dominated 
by the chirping phase, in agreement with the findings of Loren- 
Aguilar et al. (2005). Moreover, the gravitational waveforms ob- 
tained here are very similar to those computed by Loren-Aguilar 
et al. (2005) and, thus, do not depend appreciably on the num- 
ber of particles used to calculate them. This stems from the fact 
that most of the emission of gravitational waves comes from 
the regions in which the SPH particles change appreciably their 
velocities, and these regions were well resolved in both sets of 
simulations. Larger order terms of gravitational wave emission 
could be included in the calculation of the strains. These terms 
include the current-quadrupole and the mass octupole. It has 
been shown (Schutz & Ricci 2001) that for the first of these to 
be relevant an oscillating angular momentum distribution with 
a dipole moment along the angular momentum axis is needed. 
Consequently, in our calculations only the mass octupole should 
be considered in the best of the cases. Within this approximation 
a term of the order of v/c ~ 10"^ would be added to the derived 
strains. We have performed a post-processing of our simulations 
and we have found that the octupole emission is of the order of 
h ~ 10"^'*, just above the numerical noise (h ~ 10"^^), but in 
any case totally negligible. Since the gravitational wave signal 
is dominated by that of the inspiralling phase, in order to assess 
the feasibility of detecting it using gravitational wave detectors 
we have assumed that the orbital separation of the binary sys- 
tem is that of the last stable orbit. Furthermore, we have also 
assumed that the integration time of LISA will be one year It 
is then straightforward to demonstrate that during this time in- 
terval the variation of the orbital separation is negligible. With 
these assumptions the double white dwarf binary system basi- 
cally radiates a monochromatic wave and it is easy to assess the 
feasibility of detecting the signal produced by the coalescence 
of close binary white dwarf systems. This is done in figure |7] 
where we show the strength of the resulting signals and we com- 
pare them with the spectral distribution of noise of LISA, when 
a distance of 10 kpc is adopted. As can be seen in this figure all 
the systems are well inside the detectability region and, hence, 
LISA should be able to distinguish them from Galactic noise. 
However, the typical rate of white dwarf mergers is rather small, 
of the order of ~ 9 x 10"^ yr ' (Nelemans 2003) and, conse- 
quently, although there is an uncertainty of a factor of 5 in the 
rate of white dwarf mergers, the expected detection rate is small 
as well. 



Qj^(t-R) ^ Piju(N)J]m(p)[2v\p)v'(p) 



+ x''ip)a'ip) + x'(p)a'ip)] 



where 



Pijki(N) = {6ij - NiNkXdj, - NjN,) 
- ]:iSij-NiNjmi-NkN,) 



(15) 



(16) 



is the transverse-traceless projection operator onto the plane or- 
thogonal to the outgoing wave direction, A^, m(p) is the mass of 
each SPH particle, and x{p), v(p) and a(p) are, respectively, its 
position, velocity and acceleration. 



5.3. Fallback luminosities 

Another potential observational signature of the mergers studied 
here is the emission of high-energy photons from the fallback 
material in the aftermath of the coalescence itself. We have al- 
ready shown that as a result of the merger of two white dwarfs 
of different masses, most of the SPH particles of the disrupted 
secondary form a keplerian disk. These SPH particles have cir- 
cularized orbits. However, as it occurs in the coalescence of dou- 
ble neutron stars, some material of the secondary is found to be 
in highly eccentric orbits as well. After some time, this mate- 
rial will most likely interact with the recently formed disk. As 
discussed in Rosswog (2007) the timescale for this is not set by 
viscous dissipation but, instead, by the distribution of eccentri- 
ties. We follow closely the model proposed by Rosswog (2007) 
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dwarf binary systems studied here when, a distance of 10 kpc 
is adopted, with the spectral distribution of noise of LISA. The 
spectral distribution of noise of LISA is for a one year integration 
period. We have adopted a signal-to-noise ratio rj = 5. 
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Fig. 8. Fallback accretion luminosity for our three fiducial dou- 
ble white dwarf mergers. The units of time are seconds, whereas 
those of luminosities are erg/s. A straight line with slope 5/3 is 
shown for the sake of comparison. 



and calculate the accretion luminosity obtained from the inter- 
action of the stellar material with high eccentricities with the 
newly formed disk by assuming that the kinetic energy of these 
particles is dissipated within the radius of the debris disk. 



In figure [8] we have plotted the accretion luminosities as a 
function of time for our three fiducial cases. We emphasize that 
these luminosities have been computed assuming that the highly 
eccentric particles loose all its kinetic energy when interacting 
with the disk, for which we adopt the radius obtained by the 
end of our SPH simulations, which are those shown in Table [T] 
Moreover, only a fraction of this energy will be released in the 
form of high-energy photons. Thus, the results shown in Fig. 
[8] can be regarded as an upper limit for the actual luminosity 
of high-energy photons. Note that although the luminosities are 
smaller than those typically obtained for the merger of double 
neutron stars — which are typically of the order of ~ 10^^ erg/s 
— white dwarfs mergers predict a very similar time dependence 

This is an important result because it shows that observa- 
tions of high-energy photons can help in detecting the gravita- 
tional wave signal radiated by these systems. In fact, the detec- 
tion of the gravitational waves arising from the merger of white 
dwarfs is a tough task because, as previously explained, the sig- 
nal is dominated by the inspiralling phase and the waveforms do 
not have a prominent peak before the ringdown phase. Thus a 
combined strategy in which optical, UV, X-ray and gravitational 
wave detectors are used could be very useful. 

5.4. Long-term evolution 

We have already shown that non-explosive nuclear burning takes 
place during the merging phase. However, this does not necessar- 
ily mean that such an explosion could not take place due to mass 
accretion from the disk at late times. If mass acccretion occurs at 
rates smaller than 10"^ Mq yr"' then, central carbon ignition is 
possible and a SNIa is the most probable outcome. On the other 
hand, if the accretion rates are larger than this value, then off- 
center carbon ignition is the most probable outcome, giving rise 
to an inward propagating burning flame and an ONe white dwarf 
is likely to be formed (Nomoto & Iben 1986; Garcia-Berro 
& Iben 1994; Ritossa et al. 1999) wich might eventually form 
a neutron star by accretion-induced collapse (Saio & Nomoto 
1985; Woosley & Weaver 1986; Gutierrez et al. 1996, 2005). 
However, once the disk has been formed, angular momentum 
viscous transfer is relevant and the hydrodynamical timescale 
of the disk becomes very large. Consequenlty, the subsequent 
evolution of the disk cannot be followed using a SPH code. 
However, some estimates of the accretion rate can still be done 
by considering the typical viscous transport timescales. 

The typical viscous transport timescale is (Mochkovitch & 
Livio 1989, 1990) 



1 dT 
f~dt 



where T is the rotational kinetic energy and 



dt 



(17) 



(18) 



is its rate of change. In this expression Q is the angular veloc- 
ity, r is the radial cylindrical coordinate and rj is the (physical) 
viscosity parameter, that depends on the viscous mechanism. If 
the disk is laminar and the viscosity is that of degenerate elec- 
trons then rj - 2.0 x 10"^p^^^g/cm s (Durisen 1973, Itoh et al. 
1987) and the associated accretion rates can be obtained taking 
into account that 



GM,M 



(19) 
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Fig. 9. Richardson number as a function of the distance. When 
Ri> 1/4 the disk is turbulent. The horizontal dashed line corre- 
sponds to Ri=l/4, 



where M» and /?» are the mass and radius of the central object. 
If, instead, the disk is turbulent the classical approximation of 
Shakura & Sunyaev (1973) is valid. Within this approximation 
the viscous timescale is given by 



-1 /^disk W^disk 

I H I 



(20) 



where a ~ 0.1 is the standard viscosity, Cj is the sound speed, 
Rdisk is the radius of the disk and H is the disk half-thickness. 
Both the radius of the disk and its half-thickness are listed in 
Table [T] for each one of the simulations presented here. The ac- 
cretion rate is then given by 



M =^ 



Mdisk 



(21) 



In order to check if the disk is turbulent we have computed 
the Richardson number 



Ri = 



(g)(l-f) 



(22) 



In this expression g^^ is the effective gravity, that is the real 
gravity minus the centrifugal force, y is the logarithmic deriva- 
tive of the pressure with respect to the density, F is the adiabatic 
index and the rest of the symbols have their usual meaning. We 
have chosen y = 1.4 and F = 5/3. If the Richardson number 
is smaller than 1/4 stability against turbulence is guaranteed. If 
this is not the case the disk may be turbulent, since this is only a 
necessary condition, but Ri> 1 /4 is a good indication for turbu- 
lence to occur In figure |9] we show the Richardson number as a 
function of the radial coordinate for the three fiducial cases de- 
scribed here. As can be seen in this figure the condition Ri> 1/4 
is satisfied in the innermost regions of the disk, up to distances 
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Table 3. Typical viscous timescales (in years), dissipated kinetic 
energies (in erg/s), and accretion rates (in Mq yr"') for the case 
in which the laminar viscosity is used. Columns 5 and 6 show the 
corresponding viscous timescales and accretion rates obtained 
when the classical Shakura & Sunyaev (1973) expression is em- 
ployed. 



~ 0.2 Rq. Thus, the innermost regions of the disk are potentially 
turbulent and it is likely that the accretion rate is that given by 
Eq. (21 1. Nevertheless, Eq. ( 17 i provides a safe lower limit for 
the typical transport timescale of the disk. The total rotational ki- 
netic luminosity dissipated, t, and the corresponding timescales 
using this approach are given in the left section of table[3] For the 
case in which the laminar viscosity is used the resulting accretion 
rates turn out to be $ lO"'^M0/yr, and consequently central 
carbon ignition leading to a SNIa is possible. When the classical 
Shakura & Sunyaev (1973) expression is adopted the accretion 
rates are shown in the right section of table [3] As can be seen 
these accretion rates turn out to be very large. There are experi- 
mental and theoretical reasons to suspect that the central object 
will not be able to accrete material from the surrounding disk at 
these very high accretion rates. From a theoretical perspective 
it is rather clear that these accretion rates exceed the Eddington 
limit, which is order of lO^^MQyr"'. Additionally, and from 
a experimental point of view, there is growing evidence (Ji et 
al. 2006) that hydrodynamic turbulence cannot transport angular 
momentum effectively in astrophysical disks, even at very large 
Reynolds numbers, leaving as the only possible way to lose an- 
gular momentum the magnetorotational instability. 

Yoon et al. (2007) have systematically explored the con- 
ditions for avoiding off-center carbon ignition in the merged 
configurations resulting from the coalescence of two CO white 
dwarfs. They computed the evolution of the central remnant of 
the coalescence of a 0.9 + 0.6 Mq binary white dwarf, adopting 
a realistic initial model, which includes the differentially rotat- 
ing primary, the hot envelope we also find in our simulations 
and the centrifugally supported accretion disk. Our final config- 
urations resemble very much those found by these authors and, 
consequently, the same conclusions obtained in this paper hold. 
In particular, since in our models the maximum temperature of 
the hot envelope is smaller than the carbon ignition temperature 
and the mass-accretion rate from the keplerian disk is possibly 
smaller than the critical one (M ^ 5.0x 10"* Mq yr"') it is prob- 
able that at least some of the merged configurations with a total 
mass exceeding Chandrasekhar's mass may be considered good 
candidates for Type la supernovae progenitors. 



6. Conclusions 

We have performed several high-resolution Smoothed Particle 
Hydrodynamics simulations of coalescing white dwarfs. We 
have done so for a broad range of masses and chemical composi- 
tions of the coalescing white dwarfs, which includes He, CO and 
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ONe white dwarfs. Such a parametric study using a large num- 
ber of particles (4 x 10^ SPH particles) had never been done be- 
fore. Previous works on the subject used a considerably smaller 
number — by a factor of 10 — of SPH particles (Guerrero et al. 
2004), or did not explore the full range of masses and chemical 
compositions of interest (Yoon et al. 2007) — only one merger 
was computed, that of 0.9 -i- 0.6 Mq double white dwarf. In ad- 
dition, we have included a refined treatment of the artificial vis- 
cosity. In particular, we have used an artificial viscosity formu- 
lation that is oriented at Riemann-solvers (Monaghan 1997) to- 
gether with an additional switch to suppress the excess of vis- 
cosity (Balsara 1995). With this treatment the dissipative terms 
are only applied in those regions of the fluid in which they are 
really necessary to resolve a shock. This refined treatment of the 
artificial viscosity overcomes some of the problems found in our 
previous simulations (Guerrero et al. 2004). 

In all cases, the merged configuration consists of a compact 
central object surrounded by a hot corona with spheroidal shape 
and a self-gravitating keplerian disk around it. For the cases in 
which two white dwarfs of different masses are involved the re- 
sulting disk can be considered as a thin disk, whereas for the 
0.6 + 0.6 Mq case we have found that the resulting final con- 
figuration resembles a rotating ellipsoid around the central ob- 
ject with a much more modest disk. The peak temperatures at- 
tained during the merging process are smaUer than those found 
by Guerrero et al. (2004) and in line with that found by Yoon et 
al. (2007) for the case of a 0.9 + 0.6 Mq merger. We also con- 
firm the results obtained in previous works (Guerrero et al. 2004: 
Loren-Aguilar et al. 2005; Yoon et al. 2007) and we find that 
only when one of the merging white dwarfs is a He white dwarf 
nuclear reactions are relevant. However none of the cases stud- 
ied here show an explosive behavior during the merging phase. 
Furthermore, no essential differences are found when the chem- 
ical abundances obtained here using an enhanced spatial reso- 
lution and a refined prescription for the artificial viscosity and 
those obtained in previous works (Garcia-Berro et al. 2007) are 
compared. The chemical composition of the disks formed by the 
coalescence of a He white dwarf with CO white dwarfs shows an 
enrichment in heavy elements like Ca, Mg, S, Si and Fe and con- 
stitute a natural environment where planets and asteroids can be 
formed. This could explain the anomalous abundances of metals 
foimd in several hydrogen-rich white dwarfs with dusty disks 
around them and which have been attributed to the impact of as- 
teroids (Jura 2003), since it is quite unlikely that such asteroids 
could survive the red giant phase. 

We have also compared the results of our hydrodynamical 
calculations with the theoretical expectations and we have found 
a satisfactory agreement when the synchronization timescale of 
the disrupted secondary is much shorter than that of the primary. 
In this case the rate of change of the orbital distance and the cor- 
responding spins of both the donor star and of the accretor are 
reproduced with a large degree of accuracy. We have shown as 
well that the emission of gravitational waves from these kind of 
systems is strong enough to be obervable by LISA, and that the 
corresponding waveforms do not depend appreciably on the res- 
olution employed in the hydrodynamical calculations and, thus, 
that these waveforms are robust. 

We have computed as well the possible emission of high- 
energy photons produced in the aftermath of the merger. This 
high-energy emission is a consequence of the interaction of the 
material with highly eccentric orbits — which is produced dur- 
ing the first and most violent phases of the merger — with the 
resulting disk — which is formed by particles with circularized 
orbits — and we have found that the typical luminosities are of 



the order ~ 10"*' erg/s, although the precise value of the peak 
luminosity depends very much on the masses of the coalescing 
white dwarfs. The time dependence of the high-energy emission 
is oc f^/^, a behavior also present in the coalescence of double 
neutron stars and neutron stars and black holes (Rosswog 2007). 
The detection of this high-energy burst would eventually help in 
detecting the gravitational waves suposedly radiated during the 
merger. 

With respect to the long-term evolution of the mergers we 
have found that all the disks product of the coalescence of two 
white dwarfs are potentially turbulent. This result implies that 
very large accretion rates from the disk onto the primary are ex- 
pected. Despite our crude approximations, it is thus quite likely 
that these accretion rates would lead to an off-center carbon ig- 
nition, although an in depth study remains to be done. However, 
since our final configurations resemble very much those found 
by Yoon et al. (2007), it is as well probable that at least some 
of our merged configurations may be considered good Type la 
supemovae candidate progenitors. A detailed calculation of the 
evolution of the resulting disks, including an accurate descrip- 
tion of the mechanisms of angular momentum transport, must 
therefore be done. Unfortunately this task is far beyond the cur- 
rent possibilities of SPH techniques. 
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